{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/obss/BIOBSS/blob/main/examples/hrv_analysis.ipynb)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__BIOBSS - Heart Rate Variability (HRV) Analysis__\n",
    "\n",
    "_This notebook includes guidelines to help using BIOBSS for calculating heart rate variability parameters from PPG or ECG signal._"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__To run this notebook in Colab:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "git clone https://github.com/obss/biobss.git\n",
    "cd biobss\n",
    "pip install ."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!pip install -U matplotlib\n",
    "import shutil\n",
    "from biobss import FIXTURES_ROOT\n",
    "shutil.move(\"/content/biobss/sample_data\",FIXTURES_ROOT.parent)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Import required packages:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import biobss\n",
    "import os\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Table of Contents\n",
    "1. [Heart rate variability analysis using PPG signal](#hrv_ppg)<br>\n",
    "    1.1. [PPG Sample Data](#ppg_data)<br>\n",
    "    1.2. [PPG Preprocessing](#ppg_pre)<br>\n",
    "    1.3. [Calculation of HRV Parameters](#ppg_par)<br>\n",
    "2. [Heart rate variability analysis using ECG signal](#hrv_ecg)<br> \n",
    "    2.1. [ECG Sample Data](#ecg_data)<br>\n",
    "    2.2. [ECG Preprocessing](#ecg_pre)<br>\n",
    "    2.3. [Calculation of HRV Parameters](#ecg_par)<br>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Heart rate variability (HRV) analysis can be performed using either ECG or PPG signal. For this analysis, R peaks of ECG signal or systolic peaks of PPG signal should be detected first. Using peak locations, first peak to peak intervals (ppi) are derived and then HRV parameters are calculated."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "HRV parameters can be grouped into three categories as time-domain parameters, frequency-domain parameters and nonlinear parameters. BIOBSS library provides 'get_hrv_features' function to calculate HRV parameters from ECG or PPG signal. Parameters for each domain can be calculated by defining feature_types parameter or using the functions specific to each category. All functions return a dictionary of the calculated parameters.\n",
    "\n",
    "\n",
    "The HRV parameters which can be calculated using BIOBSS library are given below.\n",
    "\n",
    "<b>Time-domain HRV parameters</b>:\n",
    "- mean_nni: mean of peak to peak intervals\n",
    "- sdnn: standard deviation of peak to peak intervals\n",
    "- rmssd: root mean square of successive differences between peak to peak intervals\n",
    "- sdsd: standard deviation of successive differences between peak to peak intervals\n",
    "- nni_50: number of pairs of successive intervals that differ by more than 50 ms\n",
    "- pnni_50: ratio of nni_50 to total number of intervals\n",
    "- nni_20: number of pairs of successive intervals that differ by more than 20 ms\n",
    "- pnni_20: ratio of nni_20 to total number of intervals\n",
    "- cvnni: ratio of sdnn to mean_nni \n",
    "- cvsd: ratio of rmssd to mean_nni\n",
    "- median_nni: median of absolute values of successive differences between peak to peak intervals\n",
    "- range_nni: range of peak to peak intervals\n",
    "- mean_hr: mean heart rate\n",
    "- min_hr: minimum heart rate\n",
    "- max_hr: maximum heart rate\n",
    "- std_hr: standard deviation of heart rate\n",
    "- mad_nni: mean absolute deviation of peak to peak intervals\n",
    "- mcv_nni: ratio of mead_nni to median_nni\n",
    "- iqr_nni: interquartile range of peak to peak intervals\n",
    "\n",
    "<b>Frequency-domain HRV parameters</b>:\n",
    "- vlf: spectral power pertaining to very low frequency band (0.0033 to 0.04 Hz by default.)\n",
    "- lf: spectral power pertaining to low frequency band (0.04 to 0.15 Hz by default.)\n",
    "- hf: spectral power pertaining to high frequency band (0.15 to 0.4 Hz by default.)\n",
    "- lf_hf_ratio: ratio of lf to hf\n",
    "- total_power: sum of vlf, lf and hf\n",
    "- lfnu: normalized spectral power pertaining to low frequency band (ratio of lf to total_power)\n",
    "- hfnu: normalized spectral power pertaining to high frequency band (ratio of hf to total_power)\n",
    "- lnLF: log transformed low-frequency power\n",
    "- lnHF: log transformed high-frequency power\n",
    "- vlf_peak: max peak of power spectral density in very low frequency band\n",
    "- lf_peak: max peak of power spectral density in low frequency band\n",
    "- hf_peak: max peak of power spectral density in high frequency band\n",
    "\n",
    "<b>Nonlinear HRV parameters</b>:\n",
    "- SD1: standard deviation of Poincare plot perpendicular to the line of identity\n",
    "- SD2: standard deviation of Poincare plot along the line of identity\n",
    "- SD2_SD1: ratio of SD2 to SD1\n",
    "- CSI: cardiac stress index\n",
    "- CVI: cardiac vagal index\n",
    "- CSI_mofidied: modified cardiac stress index\n",
    "- ApEn: approximate entropy of peak to peak intervals\n",
    "- SampEn: sample entropy of peak to peak intervals"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### __Heart rate variability analysis using PPG signal__\n",
    "<a id=\"hrv_ppg\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### __PPG Sample Data__\n",
    "<a id=\"ppg_data\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "PPG sample data is provided as a csv file in BIOBSS\\sample data. The data file contains 100 PPG segments of 10-seconds length. The sampling rate is 64 Hz for all segments."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Load the sample data\n",
    "data, info = biobss.utils.load_sample_data(data_type='PPG_short')\n",
    "sig = np.asarray(data['PPG'])\n",
    "fs = info['sampling_rate']\n",
    "L = info['signal_length']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### __PPG Preprocessing__\n",
    "<a id=\"ppg_pre\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In order to detect peaks accurately, PPG signal should be filtered first. BIOBSS provides a filtering function which uses Scipy. This function can be used to implement Butterworth filter by defining the filter parameters (filter type, filter order, cutoff frequencies). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As an alternative, predefined filters can be used for each signal type. In order to use this option for PPG signal, signal_type should be selected as 'PPG'. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Filter PPG signal by using predefined filters\n",
    "filtered_ppg = biobss.preprocess.filter_signal(sig, sampling_rate=fs, signal_type='PPG', method='bandpass')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, PPG peaks should be detected to calculate interbeat (peak to peak) intervals (ibi, ppi). \n",
    "\n",
    "BIOBSS provides a peak detection function with different alternatives for the peak detection method. These methods are appropriate for PPG signal, however the parameters should be selected properly if the second peak (diastolic peak) is observable in the signal. The ___ppg_detectpeaks___ function returns a dictionary including peak locations and trough locations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Detect peaks using 'peakdet' method (delta=0.01). Delta parameter should be adjusted related to the amplitude of the signal.\n",
    "\n",
    "info=biobss.preprocess.peak_detection(sig,fs,'peakdet',delta=0.01)\n",
    "\n",
    "locs_peaks=info['Peak_locs']\n",
    "peaks=sig[locs_peaks]\n",
    "locs_onsets=info['Trough_locs']\n",
    "onsets=sig[locs_onsets]"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For further analysis, a peak control step is required to prevent errors resulting from incorrect peak detection results (missing or dupliciate peaks). The ___peak_control___ function checks the relative locations of peaks and troughs, ensuring only a single peak is located between consecutive troughs. As an option, peak correction procedure can be applied by setting the parameter __correct_peaks__ of ___ppg_detectpeaks___ function as True."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Correct the peak detection results by considering the order of peaks and troughs\n",
    "info=biobss.ppgtools.peak_control(sig=sig, peaks_locs=locs_peaks, troughs_locs=locs_onsets)\n",
    "\n",
    "locs_peaks=info['Peak_locs']\n",
    "peaks=sig[locs_peaks]\n",
    "locs_onsets=info['Trough_locs']\n",
    "onsets=sig[locs_onsets]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### __Calculation of HRV Parameters__\n",
    "<a id=\"ppg_par\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "BIOBSS provides ___get_hrv_features___ function to calculate HRV parameters from PPG signal. This function can take 'ppi', 'peaks' or 'troughs' as argument for input type. If 'ppi' is selected as input_type, ppi intervals must be provided in milliseconds. If 'peaks' or 'troughs' is selected as input_type, ppi is calculated by the function itself using peak/trough locations and sampling rate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate HRV parameters using ppi intervals\n",
    "ppi = 1000*np.diff(locs_peaks)/fs \n",
    "ppg_hrv = biobss.hrvtools.get_hrv_features(sampling_rate=fs, signal_length=L, input_type='ppi',ppi=ppi)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate HRV parameters using peak locations\n",
    "ppg_hrv = biobss.hrvtools.get_hrv_features(sampling_rate=fs, signal_length=L, input_type='peaks',peaks_locs=locs_peaks)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate HRV parameters using trough locations\n",
    "ppg_hrv = biobss.hrvtools.get_hrv_features(sampling_rate=fs, signal_length=L, input_type='troughs',troughs_locs=locs_onsets)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The functions ___hrv_time_features___, ___hrv_freq_features___, ___hrv_nl_features___ can be used seperately in order to calculate hrv parameters for a specific domain. Note that, in order to use these functions peak to peak intervals should be provided as input in milliseconds."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate time-domain HRV parameters \n",
    "ppg_hrv_time = biobss.hrvtools.hrv_time_features(ppi=ppi, sampling_rate=fs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate frequency-domain HRV parameters \n",
    "ppg_hrv_freq = biobss.hrvtools.hrv_freq_features(ppi=ppi, sampling_rate=fs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate nonlinear HRV parameters \n",
    "ppg_hrv_nl = biobss.hrvtools.hrv_nl_features(ppi=ppi, sampling_rate=fs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### __Heart rate variability analysis using ECG signal__\n",
    "<a id=\"hrv_ecg\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### __ECG Sample Data__\n",
    "<a id=\"ecg_data\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ECG sample data is provided as a csv file in BIOBSS\\sample data. The data file contains an ECG signal of 5 minutes length, sampled at 256 Hz."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Load the sample data\n",
    "data, info = biobss.utils.load_sample_data(data_type='ECG')\n",
    "sig = np.asarray(data['ECG'])\n",
    "fs = info['sampling_rate']\n",
    "L = info['signal_length']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### __ECG Preprocessing__\n",
    "<a id=\"ecg_pre\"></a>"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In order to detect R peaks accurately, ECG signal should be filtered first. BIOBSS provides a filtering function which uses Scipy. This function can be used to implement Butterworth filter by defining the filter parameters (filter type, filter order, cutoff frequencies)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As an alternative, predefined filters can be used for each signal type. In order to use this option for ECG signal, signal_type should be selected as 'ECG'. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Filter ECG signal by using predefined filters\n",
    "filtered_ecg=biobss.preprocess.filter_signal(sig, sampling_rate=fs, signal_type='ECG', method='pantompkins')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, ECG R peaks should be detected to calculate peak to peak intervals (ppi, rri). "
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "BIOBSS provides a specialized peak detection function for ECG signal. The function uses __ecgdetectors__ package and returns a dictionary of R-peak locations and amplitudes. The available methods are: 'pantompkins', 'hamilton' and 'elgendi'. \n",
    "For more, see: https://github.com/berndporr/py-ecg-detectors/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Detect R-peaks of ECG signal using 'pantompkins' method.\n",
    "locs_peaks=biobss.ecgtools.ecg_detectpeaks(filtered_ecg,fs,'pantompkins')\n",
    "peaks = filtered_ecg[locs_peaks]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### __Calculation of HRV Parameters__\n",
    "<a id=\"ecg_par\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "BIOBSS provides ___get_hrv_features___ function to calculate HRV parameters from ECG signal. When the signal type is selected as 'ECG', the function can take 'ppi' or 'peaks' as argument. If 'ppi' is selected as input_type, ppi intervals must be provided in milliseconds. If 'peaks' is selected as input_type, ppi is calculated by the function itself using peak locations and sampling rate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate HRV parameters using peak to peak intervals\n",
    "rri = 1000*np.diff(locs_peaks)/fs\n",
    "ecg_hrv = biobss.hrvtools.get_hrv_features(sampling_rate=fs, signal_length=L, signal_type='ECG', input_type='ppi',ppi=rri)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate HRV parameters using peak locations\n",
    "ecg_hrv = biobss.hrvtools.get_hrv_features(sampling_rate=fs, signal_length=L, signal_type='ECG', input_type='peaks',peaks_locs=locs_peaks)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The functions ___hrv_time_features___, ___hrv_freq_features___, ___hrv_nl_features___ can be used seperately in order to calculate hrv parameters for a specific domain. Note that, in order to use these functions peak to peak intervals should be provided as input in milliseconds."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate time-domain HRV parameters \n",
    "ecg_hrv_time = biobss.hrvtools.hrv_time_features(ppi=ppi, sampling_rate=fs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate frequency-domain HRV parameters \n",
    "ecg_hrv_freq = biobss.hrvtools.hrv_freq_features(ppi=ppi, sampling_rate=fs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate nonlinear HRV parameters \n",
    "ecg_hrv_nl = biobss.hrvtools.hrv_nl_features(ppi=ppi, sampling_rate=fs)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "biolib",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.16"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "3a4bcfb23c7e6ad66c655087280fa9f4d0273121ae7909f7735d1e02563a2438"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
